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The LISA time-delay-interferometry responses to a gravitational-wave signal are rewritten in a 
form that accounts for the motion of the LISA constellation around the Sun; the responses are given 
in closed analytic forms valid for any frequency in the band accessible to LISA. We then present a 
complete procedure, based on the principle of maximum likelihood, to search for stellar-mass binary 
systems in the LISA data. We define the required optimal filters, the amplitude-maximized detection 
statistic (analogous to the T statistic used in pulsar searches with ground-based interferometers), and 
discuss the false-alarm and detection probabilities. We test the procedure in numerical simulations 
of gravitational-wave detection. 
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I. INTRODUCTION 



The Laser Interferometer Space Antenna (LISA) is a deep-space mission aimed at detecting and studying gravita- 
tional radiation in the millihertz frequency band. A joint American and European project, it is expected to be launched 
in the year 2018, and to start collecting scientific data approximately a year later, after reaching its orbital configu- 
ration of operation LISA consists of three widely separated spacecraft, flying in a triangular, almost-equilateral 
, configuration, and exchanging coherent laser beams. In contrast to ground-based, equal-arm, gravitational-wave (GW) 
interferometers, LISA will have multiple readouts, corresponding to the six laser Doppler shifts measured between 
spacecraft. Modeling each spacecraft as carrying lasers, beam splitterSjphotodetectors, and drag-free proof masses 
on each of two optical benches, Armstrong, Estabrook, and Tinto 0, [|[ 0] showed that it is possible to combine, 
with suitable time delays, the six time series of the inter-spacecraft Doppler shifts and the six time series of the 
intra-spacecraft Doppler shifts (measured between adjacent optical benches) to cancel the otherwise overwhelming 
frequency fluctuations of the lasers (Av/u ~ 10~ 13 /VHz), and the noise due to the mechanical vibrations of the 
optical benches (which could be as large as Avjv ~ 10~ 16 /vTIz). The strain sensitivity level that then becomes 
q-i achievable, h ~ 10 _21 /\/Hz, is set by the buffeting of the drag-free proof masses inside each optical bench, and by 
the shot noise at the photodetectors. Several such laser-noise-free interferometric combinations are possible, and they 
show different couplings to gravitational waves and to the remaining system noises @, 0, 0, H| • The technique used to 
synthesize these combinations is known as time-delay i9nterfcromctry (TDI) ; in the case of a stationary array .it was 
shown that the space of all the possible TDI observables can be constructed by combining four generators @, Q ■ 

Recently, it was pointed out 0, d, 0, [H| that the rotational motion of the LISA array around the Sun and the 
time dependence of light travel times introduced by the relative (shearing) motion of the spacecraft have the effect 
of preventing the suppression of laser frequency fluctuations, at least under the current stability requirements, to the 
level of the secondary noises in the TDI observables as derived for a stationary array. This problem was addressed 
by devising new combinations that are capable of suppressing the laser frequency fluctuations below the secondary 
noises for a rotating LISA array @, [1], and for a rotating and shearing LISA array d, [Tof . In this context, the 
original stationary-array combinations are sometimes known as "TDI 1.0" (or first- generation TDI), the rotating- 
LISA combinations as "TDI 1.5," and the rotating and shearing-LISA combinations as "TDI 2.0;" following Ref. 
[To| . we refer to the last as second-generation TDI. Second-generation combinations are essentially finite differences of 
first-generation combinations, and as such they appear more complicated. However, they retain the same sensitivity 
to incoming GWs: this is because the corrections introduced in the original combinations by the changing array 
geometry are obviously important for laser frequency fluctuations, but they are negligibly small for the GW responses 
and for the secondary noises; thus, once laser frequency noise is removed, the second-generation observables become 
finite differences of the corresponding first-generation observables. At a fixed frequency, the ratio of GW response to 
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secondary noises (and hence the sensitivity) is then unchanged. 

The GW responses of the TDI combinations depend on the relative orientation of the LISA array with respect to 
the direction of propagation of the GW signal, on the strength and polarization of the signal, and on its frequency 
components. Analytic expressions for the TDI responses were first derived by Armstrong, Estabrook and Tinto [2|, 
for a stationary LISA array. A realistic model of LISA must however include the motion of the array around the 
Sun, which introduces slow modulations in the phase and amplitude of the GW responses (in addition, of course, to 
the modifications introduced by adopting second-generation TDI). For instance, the LISA responses to the sinusoidal 
signal emitted by a binary system are not simple sinusoids, but rather superpositions of many sinusoids of smaller 
amplitude. To maximize the likelihood of source detection, these effects must be modeled in GW search algorithms, 
either by including the modulations in the theoretical models of the signals (i.e., the templates) , or by demodulating 
the LISA data for a given set of sky positions as the first step of data analysis [n], [l2| • 

In this paper we derive the response of the second-generation TDI observables to the GW signals generated by a 
binary system, and we describe how signal templates based on these responses can be used in a maximum-likelihood, 
matched-filtering framework to search for binaries and to estimate their parameters once they are found. Other 
methods to analyze the LISA data for signals from binaries, implemented in the long-wavelength approximation, 
have been proposed in Refs. [H, We work in the solar-system-baricentric frame, and we follow closely the 

derivation given by Jaranowski, Krolak, and Schutz [l5j for continuous sources and ground-based detectors. A similar 
formalism was used by Giampieri [l6j to obtain the antenna pattern of an arbitrary orbiting interferometer, in the 
long-wavelength approximation. The response of an orbiting equal-arm-Michelson interferometer to a sinusoidal 
signal was worked out by Cutler 17 1, again in the long- wavelength limit. Seto [l8| extended Cutler's formalism to 
high frequencies (and to noise-canceling observables), in the context of studying optimal- filtering parameter estimation 
for supermassive-black-hole binaries. Cornish, Rubbo, and Poujade [H^, |20| obtained general expressions valid in the 
entire LISA frequency band, and for arbitrary GW signals; these expressions are given as integrals over the LISA 
arms, and they provide the basic building blocks to assemble the TDI observables. By contrast, in this paper we 
work out explicit time-domain expressions for the LISA response to moderately chirping binary systems, for all the 
second-generation TDI combinations. These expressions arc valid over the entire LISA frequency band, and they are 
written as linear combinations of four time-dependent functions; this linear structure facilitates the computation of 
matched filters and the design of optimal filtering algorithms. 

This paper is organized as follows. In Sec. II we give a brief overview of the derivation of the TDI responses to 
GWs for a stationary array, and we argue that the corrections introduced by the motion of the LISA array and by the 
time dependence of light travel times are negligibly small. Working in the solar-system-baricentric frame, we obtain 
general expressions for the GW responses of the Michelson (Xi, X 2 , A 3 ), Sagnac (ai, a 2 , a 3 ), and optimal (A, E, T; 
see [2l|) second-generation TDI observables (expressions for the first-generation observables are given in Appendix 
[Cjk finally, we derive the corresponding closed- form analytic expressions for moderately chirping binary systems, valid 
at any frequency in the LISA band. In Sec. Ill we provide expressions for the spectral densities of noise in the TDI 
combinations. In Sec. IV we combine the results of Sees. II and III to design optimal filters that can be applied to 
the LISA TDI data to search for binary stars; we take advantage of the linear structure of the responses to define an 
optimal detection statistic that does not depend on the effective polarization and on the initial phase of the binary, in 
analogy to the T statistic [TH, [22[ used in searches for continuous GW sources with ground-based interferometers. In 
Sec. V we derive the false-alarm and false-dismissal probabilities for our LISA T statistic. Last, in Sec. VI we describe 
an efficient algorithm to compute J 7 , and we implement it numerically; we perform a simulation of GW detections in 
both the low and high-frequency part of the LISA band, and for both Michelson and optimal [2l[ TDI observables, 
and we show that our algorithm yields very accurate estimates of source parameters. In the rest of this paper we shall 
use units where c = 1. 



II. TIME-DELAY INTERFEROMETRY 



Figure [T] shows the overall LISA geometry. The spacecraft are labeled 1, 2, and 3; the arms are labeled with the 
index of the opposite spacecraft (e.g., arm 1 lies between spacecraft 2 and 3). The light travel time (or, loosely, 
the armlength) along arm i is denoted by L t 0, H, H, E^l- 

The basic constituents of the TDI observables are the 
time series of the relative laser- frequency fluctuations measured between spacecraft, which are denoted by yij(t), with 
i =/= j: for instance, 2/31 (i) is the time series of relative frequency fluctuations measured for reception at spacecraft 
1 with transmission from spacecraft 3 (along arm 2); similarly, y 2 i(t) is the time series measured for reception at 
spacecraft 1 with transmission from spacecraft 2 (along arm 3), and so on. Six more time series result from comparing 
the laser beams exchanged between adjacent optical benches within each spacecraft; these time series are denoted by 
Zij, with i,j = 1,2,3, i / j (see 0, 0, E3] for details). Delayed time series are denoted by commas: for instance, 
2/31,2 = 2/31 - L 2 (t)}, and so on. 
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FIG. 1: Schematic LISA configuration. The spacecraft are labeled 1, 2, and 3; each spacecraft contains two optical benches, 
denoted by 1, 1*, as indicated. The optical paths are denoted by Li, where the index i corresponds to the opposite 
spacecraft. The unit vectors &i point between pairs of spacecraft, with the orientation indicated. 

The frequency fluctuations introduced by the lasers, by the optical benches, by the proof masses, by the fiber optics, 
and by the measurement itself at the photo-detector (i.e., the shot-noise fluctuations) enter the Doppler observables 
yij and Zij with specific time signatures; see Refs. 0,111, HH for a detailed discussion. The contribution y^ w due to 
GW signals was derived in Ref. 0] in the case of a stationary array. (Note that in Ref. 0], and indeed in all the 
literature on first-generation TDI, the notation indicates the one-way Doppler measurement for the laser beam 
received at spacecraft j and traveling along arm i. In this paper we conform to the notation used in Refs. 7, 8, 9, 10]). 

Since the motion of the LISA array around the Sun introduces a difference between (and a time dependence in) the 
corotating and counterrotating light travel times, the exact expressions for the GW contributions to the various first- 
generation TDI combinations will in principle differ from the expressions valid for a stationary array 0]. However, 
the magnitude of the corrections introduced by the motion of the array are proportional to the product between 
the time derivative of the GW amplitude and the difference between the actual light travel times and those valid 
for a stationary array. At 1 Hz, for instance, the larger correction to the signal (due to the difference between the 
corotating and counterrotating light travel times) is two orders of magnitude smaller than the main signal. Since the 
amplitude of this correction scales linearly with the Fourier frequency, we can completely disregard this effect (and the 
weaker effect due to the time dependence of the light travel times) over the entire LISA band [l(|. Furthermore, since 
along the LISA orbit the three armlengths will differ at most by ~ l%-2%, the degradation in signal-to-noise ratio 
introduced by adopting signal templates that neglect the inequality of the armlengths will be at most a few percent. 
For these reasons, in what follows we shall derive the GW responses of various second-generation TDI observables 
by disregarding the differences in the delay times experienced by light propagating clockwise and counterclockwise, 
and by assuming the three LISA armlengths to be constant and equal to L = 5 x 10 6 km ~ 16.67s [23T ]. These 
approximations, together with the treatment of the moving-LISA GW response discussed at the end of Sec. Ill C| are 
essentially equivalent to the rigid adiabatic approximation of Ref. [2Q] , and to the formalism of Ref. [HI . 



We denote the positions of the three spacecrafts by p, and the unit vectors along the arms by n^, where fii points 
from spacecraft 3 to 2, ri2 points from spacecraft 1 to 3, and 113 points from spacecraft 2 to 1. In the coordinate frame 
where the spacecraft are at rest, we can set without loss of generality 



A. Geometry of the orbiting LISA array 



pf = (L/v / 3)( 



cos 2o~i, sin 2<7i, 0), 



(1) 



and 



hf = (cosCTj, sincri, 0) 



(2) 



where 



er 4 = 3tt/2- 2(i - 1)tt/3. 



(3) 
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FIG. 2: Orbital motion of the LISA detector, shown in a solar-system-baricentric ecliptic coordinate system. The trajectories 
shown correspond to setting £ = — 7r/6, Q = 2n/yT, and £o = ijo = in Eqs. (J5J) and (|SJ. 



Because the motion of the LISA guiding center (i.e., the baricenter of the formation) is contained in the plane of the 
ecliptic, it is convenient to work in a solar-system-baricentric (SSB) ecliptic coordinate system. We take the x axis 
of this system to be directed toward the vernal point. A realistic set of orbits for the spacecraft (23|, shown in Fig. [2l 
is obtained by setting 

Pi (t)=v(t)+0 2 -pf, n i (t) = 2 -nf, (4) 

where r is the vector from the origin of the SSB coordinate system to the LISA guiding center, as described by the 
SSB components 

r = R(cosr), sin 77, 0), R = 1 AU; (5) 

the function 7/ = Qt + 770 returns the true anomaly of the motion of the LISA guiding center around the Sun. The 
rotation matrix O2 models the cartwheeling motions of the spacecraft along their inclined orbits, shown in Fig. [3j it 
is given by 

sin 77 cos £ — cos 77 sin £ sin £ — sin 77 sin £ — cos 77 sin £ cos £ — cos 77 cos £ 
02= ( — cos 77 cos £ — sin 77 sin £ sin £ cos 77 sin £ — sin 77 sin £ cos £ —sin 77 cos £ | ; (6) 
cos £ sin £ cos £ cos £ — sin£ 

the function £ = — fit + £ n returns the phase of the motion of each spacecraft around the guiding center, while £ sets 
the inclination of the orbital plane with respect to the ecliptic. For the LISA trajectory, = 27r/yr and £ = — it/6 
[2^ . For simplicity, we can set 770 = £0 = 0, so that at time t = the LISA guiding center lies on the positive x axis of 
the SSB system, while pi lies on the negative y axis. The spacecraft orbits described by Eq. |4} can be approximately 
mapped to those used by Cornish and Rubbo L 19l| by identifying our spacecraft 1,2, and 3 with their spacecraft 0, 2, 
and 1, and by setting 770 = K, £0 = 37r/2 — k + A, where k and A are the parameters defined below Eqs. (56) and (57) 
of Ref. 0J. 



B. Generic plane waveform 

At the origin of the SSB frame, the transverse-traceless metric perturbation due to a source located at ecliptic 
latitude (3 and longitude A can be written as 

H(t) = d • H s (t) ■ 0^, (7) 
where the metric perturbation in the source frame is taken to be 

h+(t) h x (t) 

H i '(t)= ( -h+(t) ) , (8) 
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FIG. 3: Cartwheeling motion of the LISA array, as plotted in a frame with center in the LISA guiding center and axes parallel 
to the SSB ecliptic frame. We show three snapshots at different times along the LISA orbital period, I yr. 




FIG. 4: Conventional definition of the GW polarizations + (dashed) and x (solid) for various ecliptic latitudes f3 and longitudes 
A. 



with h + (t) and h x (t) the two GW polarizations, and where 

(sin A cos ip — cos A sin (3 sin ip — sin A sin ip — cos A sin (3 cos ip 
— cos A cos ip — sin A sin (3 sin ip cos A sin ip — sin A sin (3 cos ip 
cos j3 sin ip cos (3 cos ip 

the dependence of the rotation matrix Oi on (3 and A enforces the transversality of the plane waves, which are 
propagating from a source located in the direction 

k = (cos A cos (3, sin A cos (3, sin f3); (10) 

the polarization angle ip encodes a rotation around the direction of wave propagation, — k, setting the convention used 
to define the two polarizations, + and x. The polarizations corresponding to ip = are shown in Fig. [4] for various 
source positions in the sky. In the center of the LISA proper frame (the frame where the spacecraft are at rest), the 
transverse-traceless metric perturbation is given by 

H L (i) = 2 - 1 (i)0 1 H s (i)Or 1 2 (t). (11) 

The time variable t that appears in h + (t) and h x (t) [and therefore in H(i) and H L (£)] is the time at the origin of the 
SSB frame. It is related to the time in the GW source frame by a relativistic time dilation, due to the proper motion 
of the source and to cosmological effects. It is however expedient to identify the two times, and to describe GW 
emission using SSB time; the time dilation is then taken into account by mapping the apparent (measured) physical 
parameters of a source into its real parameters. The source positional parameters (3, A, and ip can be mapped to the 
parameters 9, (p, and ip used in Ref. [19( by setting f3 = ir/2 — 6, A = <f>, and ip = —ip. 



— cos A cos (3 \ 

— sin A cos j3 ; (9) 
— sin j3 J 



6 



C. GW response of the LISA array 

As derived in Ref. Q for a stationary, equilateral-triangle LISA array, the one-way Doppler responses j/21 and j/31 
excited by a plane transverse-traceless GW propagating from the source direction k, are given by 

2/ 2 G i W W = [l-k.n3][* 3 (* + k-Pa-i)-*3(* + k-pi)], (12) 
yg w (t) = [l + k-A 2 ][*2(t + k.p3-L)-*2(t + k.pi)] (13) 

[in the notation of Ref. [2|, our j/21, J/31 , and k correspond to j/31, j/21, and — k respectively] where 

= - % {t \ , 2 , ^(t) = |fi;-H(i)- nj (w) 

1 - (k • rij) 2 

[the prime denotes vector transposition]. The two ^ terms in each of Eqs. (fL2"|) and (jT3]) correspond to the events of 
emission (at spacecraft 2 and 3, respectively) and reception (at spacecraft 1) of a laser photon packet; the time of the 
emission event is therefore retarded by an armlength L. The k-p^ terms represent the retardation of the gravitational 
wavefronts to the positions of the spacecraft. The other four one-way Doppler responses are obtained by cyclical 
permutation of the indices (1— » 2, 2 — > 3, 3— » 1). 

Our approximation to the GW response of the moving LISA array is obtained simply by interpreting Eqs. (|13p and 
(frjj as written in the SSB ecliptic frame, and by adopting the time-dependent equations (H]) for p^ and n,. Note that 
<&j(t) can then be written either as ®j(t) = \ h'At) ■ H(i) ■ nj(t), or $ 3 -(i) = \ (&f)' ■ H L (t) ■ hj. The time-dependent 
rotation of the fij(t) introduces an amplitude modulation of the responses, generating sidebands at frequency multiples 
of 1/yr; the time dependence of the wavefront-retardation products k • Pi(t) introduces a time-dependent Doppler 
shift caused by the relative motion of the spacecraft with respect to the SSB frame. 



D. Chirping-binary waveforms 



In the Newtonian limit, the GW signal emitted by a binary system located in the direction k can be written in the 
form of Eq. ((7|) , with 



h+{t) = /i+cos[0 s (t) + h x (t) = h* sm[<p s (t) +4,0]. 

Here (j>o is an arbitrary constant phase, and the constant amplitudes Kq and are given by 

ho = ^o(l + cos 2 t)/2, /iq = HqCOSL, 



(15) 



(16) 



where i is the angle between the normal to the orbital plane of the binary and the direction of propagation — k, and 
where 



h 



4{GM C ) 5 / 3 



c A D 



2/3 



(17) 



with Ai c = m\^m^ b / (mi+ni2) 1 ^ the chirp mass, u> the angular frequency of the GW at t = 0, and D the luminosity 
distance to the source. Last, the phase <j> s (t) is given, to the first post-Newtonian order, by [2J] 



~ - 0.(0), Mt) 



2AI 



e 5 / 8 (t) 



/3715 55 n 
V 8064 + 96M 



e- 1/4 (t) 



where 



Q(t) 



fJ,C 

5GAP 



<tc-t), 



(18) 



(19) 



tc 



GM 2 5 1 



743 
252 



924 fi s 
252 M 







~GMuj~ 


2/3 




, Xo = 


2c 3 





(20) 



with AI = mi + m,2 the total mass of the binary, and /1 
coalescence of the binary from the initial instant t = 0. 



mxmijM the reduced mass. The time t c is the time to 
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TABLE I: Contributions to the evolution of GW frequency for various types of compact, stellar-mass binaries (white dwarfs with 
m — O.35M0, neutron stars with m = 1.4Mq, and black holes with m — 6Mq), for selected (initial) GW frequencies within the 
LISA band. The contributions are expressed as GW cycles over one year of evolution, and the effects of Newtonian-order (N) 
and first post-Newtonian-order (IPN) terms are shown separately. The column labeled "Doppler" reports the integrated phase 
shift (in cycles) due to the increased Doppler shifting of the source as the frequency increases [see Eq. (|45p ], where significant. 
At / = 10" 3 Hz there is no significant evolution of GW frequency over one year. The symbol ut " indicates that the Taylor 
expansion of the phase given by Eq. (|21|l is accurate to within a quarter of a cycle. Numbers are not shown where a binary of 
a given class cannot exist at a given frequency. Some of the conclusions that can be drawn from this table are apparent also in 
Figs. 10 and 12 of Ref. [2p| ]: up to about 1 mHz, LISA cannot differentiate (using one year of data) between a monochromatic 
binary and a chirping binary (see Fig. 10 of Ref. 20]); above that frequency, chirping becomes appreciable (one additional GW 
cycle over a year in this table corresponds to a frequency shift of one bin in Fig. 12 of Ref. [I(|), but we see that it can still be 
modeled faithfully by the linear-chirp model of Eq. (I21|l . 

In Table [II for binaries consisting of various combinations of white dwarfs (WDs, with m = O.35M0), neutron stars 
(NSs, with m = 1.4M ), and black holes (BHs, with m = 6M ), and for various fiducial GW frequencies within 
the LISA band, we show the contributions to the evolution of GW frequency over one year caused by terms at the 
Newtonian (N) and first post-Newtonian (IPN) order. The table shows that at frequencies smaller or equal to 1CP 3 
Hz, the evolution of frequency is negligible. At frequencies approaching 10 mHz, the change in frequency becomes 
significant, and needs to be included in the model of the signal; however, only the first derivative of the frequency is 
needed up to about 50 mHz. In binaries with WDs of mass ~ 0.35Mq, above ~ 20 mHz the WDs fill their Roche 
lobe, and the dynamical evolution of the system is then determined by tidal interaction between the stars. In binaries 
with either a NS or a BH, post-Newtonian effects become important at about ~ 50 mHz. At 1 Hz and above, these 
binaries will coalesce in less than 1 yr; furthermore, population studies (25| suggest that the expected number of 
binaries above 50 mHz containing neutron stars and black holes is negligible. (The effects of frequency evolution in 
the LISA response to GW signals from inspiraling binaries are also discussed in Ref. [26j.) 

Therefore, for sufficiently small binary masses, for sufficiently small GW frequencies (and definitely for all non- 
tidally-interacting binaries that contain WDs), we can approximate the phase of the signal by Taylor-expanding it, 
and then neglecting terms of cubic and higher order. The resulting expression for the signal phase <f> s (t) is 

Mt)-u;t+^t 2 , where u = ^ (^j^Y ^ 11/3 . (21) 



E. TDI responses 

The response of the second-generation TDI observables to a transverse-tracele ss, p lane GW is obtained by setting 
Vij(t) = 2/ij W W [according to Eqs. (fT2"|) and fl!3[l ] in the TDI expressions of Ref. [9l.ll0|. For instance, the GW response 
of the second-generation TDI observable X\ is given by 

A l - [\VZ1 + J/13.2J + (2/21 + 2/l2,3),22 — (J/21 + 2/l2,3J ~ (.J/31 +2/13, 2>, 33j 
v * ' 

X GW (t) 

- [(2/3T + 2® + (2/2T + 2@,22 - (2/2T + 2@ " (2/3T + 2233 • (22) 

S v ' 

X GW (t-2L 2 -2L 3 )~X GW (i-4L) 

As anticipated above, here we are disregarding the effects introduced by the time dependence of light travel times, 
and by the rotation- induced difference between clockwise and counterclockwise light travel times [23] - Each of the 
two terms delimited by square brackets in Eq. (|22p corresponds to the GW response of the first-generation Michelson 
observable X 0. The TDI observables X2 and X3 are obtained by cyclical permutation of indices in Eq. (I22p . 
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Likewise, the second-generation Sagnac observables a%, a 2 , and a 3 can be written in terms of the first-generation 
Sagnac observables a, (3, and 7 (9l.[l0|: 

a GW (t) = a GW (t) - a GW (t - L x - L 2 - L 8 ) ^ a GW (t) - a GW (t - 3L). (23) 

We shall now assemble the Doppler measurements y GW from the various ingredients that enter Eqs. (fT2|) and (Tl3|) . We 
begin with the functions 3?j of Eq. (|14|) , which can be rewritten as a linear combination of the two GW polarizations 
h+(t) and h x (t): 

<t> j (t)=F;(t)h + (t) + F*(t)h»(t), (24) 

where 

Ff(t) = u j (t)cos2ip + v j (t)sm2ip, (25) 

F/(t) = Vj(t)cos2ip-Uj(t)sm2ip. (26) 

The modulation functions Ui(t) and i^(i) depend rather intricately on the LISA-to-SSB (O2) and source-to-SSB (Oi) 
rotations; thus, Ui{t) and Vi(t) depend on time through r](t) and £(i), and on the position of the source in the sky, 
given by the ecliptic coordinates /3 and A. Explicitly, we have 

Ui (t) = U cos(-2 7i ) + Ui cos(6 - 2 7i ) + U 2 cos(2(5 - 2 7i ) + U 3 cos(3<S - 2 7i ) + U 4 cos(4<S - 2 7l ) (27) 
+ (| — I cos 2 £) cos 2 /? — i sin/3 cos /3 cos ( sin £ cos 8 + \ cos 2 C(l — \ cos 2 /?) cos 25, 

Vl (t) = V sin(-2 7i ) + Vi sm(S - 2 7i ) + V 3 sin(3<S - 2 7l ) + V 4 sin(4£ - 2 7l ) (28) 
— ^ cos £ sin £ cos /3 sin 6 + \ cos £ sin (3 sin 2£, 

where 

5(t) = A - ?7(t) = A - 770 - fit, (29) 
7i = A - 770 - Co ~ (30) 

[see Eq. ([3]) for the definition of Ci, and remember that rj = £lt + i]o, £, = — fit + ^q], and where the coefficients Ui and 
Vj are given by 

U = ^(l + sin 2 /3)(l-sinC) 2 , (31) 

Ui = -|sin2/?cosC(l-sinC), (32) 

U 2 = § cos 2 cos 2 C, (33) 

[7 3 = ism2/3cosC(l + smC), (34) 

Ui = ^(l + sin 2 /3)(l + sinC) 2 , (35) 

V = -|sin/3(l -sinC) 2 , (36) 

Vi = |cos/3cosC(l -sinC), (37) 

V 3 = icos/3cosC(l + sinC), (38) 

Vi = |sin/3(l + sinC) 2 , (39) 

with £ = — 7r/6. Expanding the antenna patterns F^(t) and F*(t) of Eq. (|24| . and using trigonometric identities to 
absorb the initial phase 0o into constant coefficients, the functions &j(t) can be finally written as 

= a (1) u.,(i)cos0 s (i) + a (2) v 3 (t) cos (f> s (t) + a (3) ttj(i) sin0 s (t) +o (4) j) i (t)sin^(t), (40) 

where the constant amplitudes a( fc ) are given by 



a 



(1) 



/ij cos 0o cos 2-0 — Zip sm <^o sm 2V 1 , (41) 

a*- 2 -* = cos 0o sin 2ip + sin 0o cos 2ip, (42) 

a*- 3 -* = — /ijsin0o cos 2-0 — Hq cos 0o sin 2-0, (43) 

a*- 4 - 1 = — /ij sin 0o sin 2-0 + Zip cos 0o cos 2-0. (44) 
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Because the timescale of detector motion is much longer than the typical GW period (and because we are neglecting 
the evolution of the GW amplitude ho), it is sufficient to apply the retardations k • Pi(t) of Eqs. (fT!2)) and (fT5)) to the 
GW phase: 

<f> s (t + k • ~ ut + \iut 2 + (uj + tut) k ■ [r + 2 (i) ■ pf ] 

~ tut + \ujt 2 + (u + Cut) i?cos/3cos(Oi + rjo - A) +tuk- (0 2 (i) ■ pf) (45) 
* , ' v . ' 

<f>(t) 2uLdi(t) 

where <j>(t) is the GW phase retarded to the position of the LISA guiding center, and where we have defined 

di(t) = k ' 'Pf) = ^ cos/3sin724 _ i s in^siii(<J-7 2i ) + ^|cos/3sin(2J-72i) (46) 

(setting £ = —ir/6, and with 72 ^ = A — 770 — £0 + 2(7,;). Equations (fT2")l and (fT5|) contain also the projection factors 
k • n; (t) , which are given explicitly by 

d(t) = -k- [0 2 (t) -nf] = |cos/3sin7i - ^ sin/?sin(<S - 7l ) + ± cos/3sin(25 - 7i ). (47) 

The functions Cj(t) and di(t) are related by d\ = (c 2 — c 3 )/6, d 2 = (C3 — ci)/6, and <i 3 = (ci — c 2 )/6. Substituting the 
expressions (fT5)) [and similar ones] for the y§ w into Eq. (f2"2"| for Xp w , we get after some algebra 

4 

Xf w (i) = 4wL wa.(uL) sin(2wL)£V*>x{ fc) (t); (48) 

fe=i 

the functions (i) are given by 

{sinc[(l + c 2 )x/2] sin [<£(*) - xd 2 - 7x/2] + sinc[(l - c 2 )x/2] sin [</>(*) - xd 2 - 9x/2] } 

{sinc[(l + c 3 )x/2] sin[0(i) - a;d 3 - 9x/2] +sinc[(l - c 3 )x/2] sin - a;d 3 - 7x/2]},(49) 

{sinc[(l + c 2 )x/2] cos[^(t) - xd 2 - 7x/2] + sinc[(l - c 2 )x/2] cos[0(*) - zefe - 9x/2] } 

{sinc[(l + c 3 )x/2] cos[0(i) - xd 3 - 9x/2] + sinc[(l - c 3 )x/2] cos[0(i) - xd 3 - 7x/2] }{50) 

where x — loL and sine . . . = (sin . . .)/(. . .). The GW responses for X 2 and X3 can be obtained by cyclical permutation 
of the spacecraft indices. 

The GW response for a% can be written in similar form: 

fe=4 

af w (t) = 2^Lsin(f W L)5>«c4 fc) (i), (51) 

fe=l 







«a(t) 


Jffl 




«2(*) 



«a(t) 
«3(*) 







Ua(t) 









U3(*) 
«3(t) 
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where 







Ul(i) 


(2) 




Vi(t)_ 



u 2 (t) 
«a(*) 

«s(*) 

«3(*) 







Ul(f) 


(4) 




_«l(fj_ 



ua(t) 
«s(*) 

«3(*) 



{sinc[(l + ci)x/2] cos [0(f) — xdi — 3x] — sinc[(l — ci)x/2] cos [0(f) — xdi — 3x]} 

{sinc[(l + c 2 )x/2] cos [0(f) - xd 2 - 2x] - sinc[(l - c 2 )x/2] cos [0(f) - xd 2 - 4x] } (52) 
{sinc[(l + 03)35/2] cos [0(f) — xd 3 — Ax] — sinc[(l — 03)^/2] cos [0(f) — xd 3 — 2x] }, 
{sinc[(l + c\)x/2\ sin [0(f) — xdi — 3x] — sinc[(l — ci)x/2] sin [0(f) — xdi — 3x] } 

{sinc[(l + C2)af/2] sin[0(f) - xd 2 - 2x] - sinc[(l - c 2 )x/2] sin[0(f) - xd 2 - Ax] } (53) 
{sinc[(l + 03)35/2] sin [0(f) — xd% — Ax] — sinc[(l — 03)35/2] sin [0(f) — xd 3 — 2x] }. 



The a 2 and 013 combinations are again obtained by cyclical permutation of the spacecraft indices. 

For the second-generation TDI observable Ci (see Ref. [13]; Ci i s uniquely determined in the equal-armlength limit, 
unlike in the general case) we find 



d GW (f) = 2loL S in(i^L)^ a «Ci fe) W 



(54) 



fe=i 



with 







Ul(t) 


a(2) 
>1 




y 1 {t)_ 



u 2 {t) 

«2(*) 

u 3 (t) 

«3(*) 



>(3)" 
S>1 






a(4) 
>1 







«a(t) 
«a(*) 

u 3 (f) 

«3(*) 



{sinc[(l + ci)x/2] cos [0(f) — xdi — 3x] — sinc[(l — c\)x/2] cos [0(f) — xdi — 3x] } 

{sinc[(l + c 2 )x/2] cos [0(f) — xd 2 — 3x] — sinc[(l — c 2 )x/2] cos [0(f) — xd 2 — 3x] } (55) 
{sinc[(l + 03)^/2] cos [0(f) — xd 3 — 3x] — sinc[(l — 03)35/2] cos [0(f) — xd 3 — 3x] }, 
{sinc[(l + Ci)x/2] sin [0(f) — xdi — 3x] — sinc[(l — ci)x/2] sin [0(f) — xdi — 3x] } 

{sinc[(l + c 2 )x/2] sin [0(f) — xd 2 — 3x] — sinc[(l — 02)35/2] sin[0(f) — xd 2 — 3x] } (56) 
{sinc[(l + 03)^/2] sin [0(f) — xd 3 — 3x] — sinc[(l — C3)x/2] sin [0(f) — xd 3 — 3x] }. 



Finally, the optimal TDI observables [21| , which here we denote as A, E and T to distinguish them from the optimal 
combinations A, E, and T derived within first-generation TDI, are defined as linear combinations of ct\, a 2 and 0:3: 



A = ^("3 - Oil), 
1 . 

E = —={ai- 2a 2 + a 3 ), 
V6 

f = —^( a i + "2 + a 3 ). 



(57) 



It is clear that A, E, and T are also optimal, in the sense discussed in Ref. [2l|: this is because they can be 
written as time-delayed combinations of the first-generation optimal TDI observables, such as A = A(t) — A(t — 3L), 
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E = E{t) — E(t — 3L), and T — T(t) — T(t — 3L); since by construction the noises that enter A, E, and T are 
uncorrelated, it follows that the noises that enter A, E, and T are also uncorrelated, making these observables 
optimal. 

We recall that in the high-frequency part of the LISA band (i.e., for frequencies equal to or larger than 5 mHz), there 
exist three independent TDI GW observables (such as A, E, and T, or Xi, X 2 , and X 3 ). However, for frequencies 
smaller than 5 mHz, there are essentially only two independent observables: this is especially obvious if we reason in 
terms of the optimal combinations, where we observe that for low frequencies the GW signal response of T declines 
much faster than the responses of A and E [5j, l2l| . 



F. TDI responses in the long-wavelength limit 



The long-wavelength (LW) approximation to the GW responses is obtained by taking the leading-order terms of 
the generic expressions in the limit of 10L — > 0. For instance, for X\ [Eqs. (|^5|) — (f5U)) ] . we get 



GW 
LW 



16(wL) 3 \[u 3 {t) ~ u 2 (t)] a (1) sin0(t) - a (3) cos 0(f) + [v 3 (t) - v 2 (t)} a (2) sin (f){t) - a (4) cos <p(t) }, (58) 



with given by Eqs. (glj)-(|44|), and Ui(t), v l {t) by Eqs. {27}, (H§. The LW responses for X 2 and X 3 can be 
obtained by cyclical permutation of the indices. Adopting the notation of Ref. Q, we find also that 

X™ w (t) ~ 8L 3 [(nf )' • H L (t) • n£ - (n£)' • H L (t) ■ nf] , (59) 

where denotes the third time derivative, fif is given by Eq. ([2]), and t\ L (t) is given by Eq. (TTTj) . 

The GW responses of the Sagnac observables are equal simply to ^X^^. From Eqs. (f5T|) we then get the 

LW GW responses , E§™ , and T^. 

A£$ ~ iV2{ujLf |[2u 2 (i) - ux{t) - u 3 (t)\ a (1) sin^i) - a (3) cos0(t) 

+ [2v 2 (t) -ui(t) -v 3 (t)} a (2) sin0(i) - a {4) cos cj){t) }, (60) 



^IW 



3V6{loL) 3 \[u 3 (t) - U!(t)} a (1) sin 0(f) - a (3) cos0(t) 



[v 3 (t) a (2) sin0(t) -a (4) cos0(f) }, (61) 

(62) 



III. NOISE SPECTRAL DENSITY 



The s pec tral density of noise for the first-generation TDI observables X, Y, Z, a, /3, 7, A, E, and T is given in 
Refs. [3I l2ll| in the case of an equilateral LISA array, assuming that the noises appearing in all the proof masses and 
optical paths are uncorrelated. The finite-difference relations between first- and second-generation TDI observables 
[such as Xi(t) = X(t) — X(t — 4L), a\(t) — a(t) — a(t — 3L)] imply simple modifications to the first-generation noise 
densities: for instance, 

S Xi (lo) = 4 sin 2 (2loL)S x (lo), (63) 
S ai (lu) = 4 sin 2 (3wL/2) S a (uj); (64) 

inserting the expression of Sx = Sy = Sz from Ref. [|[ into Eq. (|6"3")) yields 

S Xl = Sx 2 = Sx 3 = 64 sin 2 (ujL) sin 2 (2uL) [2(1 + cos 2 uL)S pm + S op ] , (65) 

where S pm = 2.54 x 10~ 48 /~ 2 Hz -1 and 5 op = 1.76 x 10 -37 / 2 Hz -1 are the fractional-frequency-fluctuation spectral 
densities of proof-mass noise and optical-path noise respectively . These values correspond to an rms single-proof- 
mass acceleration noise of 3 x 10 -15 m sec -2 Hz -1 / 2 , and to an rms aggregate optical-path noise 20 x 10 -12 m Hz -1 / 2 , 
as quoted in the LISA Pre-Phase A Study [23[ . For the other TDI observables we find 

S ai = S a2 = S a3 = 8sin 2 (3wL/2)([4sin 2 (3w£/2) + 8 sin 2 (cj£/2)] S pm + 3S op ), (66) 



Sa = Se = 32sin 2 (wi/2)sin 2 (3wL/2)([6 + 4cos(wL) + 2 cos(2wi)] S* pm + [2 + cos(wL)] 5 op ), (67) 
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10" 4 10" 3 10" 2 10" 1 

Frequency [Hz] 



FIG. 5: Spectral densities of noise for the second-generation TDI observables X\ (continuous), a± (dashed), A (dotted), and T 
(dash-dotted). 

Sf = 8[l + 2cos(wi)] 2 sin 2 (3wi/2)[4sin 2 (wi/2)S' pm + S° p ] . (68) 
All the noise spectra are shown in Fig. 5. In the long-wavelength approximation, the noise expressions simplify to 

5 lw = ^lw = ^lw „ 256(uL) 2 [4(LoL) 2 S pm + {toL) 2 S op }, (69) 

S«Y = SaY = - 18(wL) 2 [ll(wL) 2 S pm + 3S op ], (70) 

s lw = 5 lw „ 5 4(uL) 2 [A(uL) 2 S pni + {u>L) 2 S op }, (71) 

S^ w ~ l62(LuL) 2 [(ujL) 2 S pnl + S op ]. (72) 



IV. OPTIMAL FILTERING OF THE LISA DATA 

In this section we develop a maximum-likelihood (ML) formalism to detect GW signals from moderately chirping 
binaries and to estimate their parameters, by analyzing the time series of the TDI observables. ML detection is based 
on maximizing the likelihood ratio A(# l ) over the source parameters 9 l ; this ratio is proportional to the probability 
that the observed detector output could have been produced by a GW source with parameters 9 l , plus instrument 
noise. The magnitude of the maximum indicates the probability that a signal is indeed present, while its location 
indicates the most likely parameters (the ML parameter estimators). Under the assumption of Gaussian, stationary, 
additive noise, log A(8 l ) is computed by correlating the detector output, x(t), with the expected GW detector response 
hgi(t), while weighting the correlation in the frequency domain by the inverse spectral density of instrument noise, 
S~ 1 (f). The family of GW responses {h S i(t)}, divided (in the frequency domain) by S n (f) to incorporate the noise 
weighting, are known as optimal filters. 

In Sec. IIV Al we describe the computation of A and of the ML parameter estimators for the optimal filters derived 
from the GW responses of Sec. II, and we show how to maximize A algebraically over the four source amplitudes a^K 
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The amplitude-maximized log A (known as J- ) is then used as a detection statistic to search for the most likely GW 
source, by maximizing it over the remaining source parameters [here denoted by thus, 9 l = (ft'*', In Sec. 
IIVBI we study the statistical distribution of in the absence (or presence) of a GW signal of parameters 

this distribution determines the statistical significance of observing a certain value of !F, for a fixed In Sec. IIV CI 
we study the statistical significance of measuring a certain value of the completely maximized statistic max^ .F(£ M ), 
which leads to the total false-alarm probability for a GW search over a range of intrinsic parameters. See Helstrom 
[28| for an extended discussion of ML detection and parameter estimation. 



A. Maximum-likelihood search method 



As discussed in Sees. I and II, the LISA Doppler measurements can be recombined into the laser-noise and optical- 
bench-noise free TDI observables, all of which can be obtained as time-delayed combinations of three generators 
Thus, in the following we denote the TDI data as the three-vector x(£) (we shall very shortly specify a convenient 
vector basis). In the case of additive noise, we can write x(i) = n(t) + h(t), where n(t) represents detector noise and 
h(t) the GW response. Idealizing n(t) as a zero-mean, Gaussian, stationary, continuous random process, we have [28[ 



logA = (x|h)-i(h|h), 



where the scalar product (. . . | . . .) is defined by 



(x|y)=4Re/ it • S" 1 • y 
Jo 



(73) 



(74) 



here the dagger denotes transposition and complex-conjugation, the tilde denotes the Fourier transform, and S denotes 
the one-sided cross spectral density matrix of detector noise, defined by the expectation value 



f[n(/)nt(/')] = |5(/-/')S„(/). 



(75) 



The larger the signal with respect to the noise, the higher the probability that a ML search (performed with 
the appropriate optimal filter) will yield a statistically significant detection, and the better the accuracy of the ML 
parameter estimators. The accuracy of estimation is also better for the parameters on which the signal is strongly 
dependent. Signal strength is characterized by the optimal signal-to-noise ratio (optimal S/N), 

P 2 = (h|h); (76) 

while the dependence of the instrument response on the parameters is characterized by the Fisher information matrix, 



ah 



dh 



(77) 



By the Cramer-Rao inequality [28| , the diagonal elements of T^- 1 provide lower bounds on the variance of any unbiased 

estimators of the l . In fact, the matrix T,^ 1 is often called the covariance matrix, because in the limit of high S/N 
the ML estimators become unbiased, and their distribution tends to a jointly Gaussian distribution with covariance 
matrix equal to r^. 1 . 

The optimal TDI observables [2l[ are obtained by diagonalizing the cross spectrum S n ; it turns out that the 
eigenvectors A, E, and T are independent of frequency. The new observables are the linear combinations of the 
Sagnac observables a\, oti, 013 given by Eq. (|57p . and by definition their noises are uncorrelated. With reference to Eq. 
(|6^|) . we define S^(f) = S^^if) = = 27r/) and Sff(f) = Sf(u) = 2irf). It is convenient to use the optimal 

observables A, E, and T as a basis for the LISA TDI observables, setting 



(78) 



the GW responses A GW , E GW , and T GW are given in Sec. Ill El for the case of moderately chirping binaries. For these 
sources, S^(u>) and Srp(uj) are approximately constant over the signal bandwidth, so we can expand Eq. (|73[) as 





A(t) 




" A GW (t) ' 


x(t) = 


E(t) 


, h(t) = 


E GW (t) 




f(t) 




f GW (t) 



log A 



T [(A\\A GW ) - i(i GW p GW ) + (E\\E GW ) 
+ T [(T||T GW ) - i(T GW ||T GW )] /Sf(u), 



±{E GW \ \E GW )] /S A (u) 



(79) 
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where To is the time of observation, and where we have introduced the time-domain scalar product 

To 



(B\\C) = (2/T ) B(t)C(t)dt. 
Jo 

Given the linear dependence [2l[ of A, E, and T on a\, a.2, and from Eq. (j5"Tj) it follows that 



(80) 



A GW (t) 

E GW (t) 
T GW (t) 



2x sin(|x) a 



(k) 



k=l 



A^(t) 



2a; sin(|a;) a 



(fc) 



fc=i 



l 



(fe) (fe; 
ai — a) 



2ai k) +a ( 3 k) 



i / (fe) (fe) , (fe) 
-t= ( a, + a\ + a' 



where once again x = wi, the amplitudes a^ fc ' are given by Eqs. (|4i p -(|44 p . and the functions a[ k ^ are given by Eqs. 
(|52p and (|53p. and by similar equations obtained by cyclical permutation of the indices. Note that the component 
functions A^(t), E^ k '(t), and T^ k \t) do not depend on the amplitudes (or equivalently, on h$ , , (f>o, and ^); 
they do however depend on the remaining (intrinsic) source parameters, u>, w, /3, and A. 

The ML parameter estimators 9 l are found by maximizing log A with respect to the source parameters 9 l : that is, 
by solving 



d log A 

86* 



0. 



(82) 



For the this is accomplished easily by solving the linear system 



^M««a«=iV« 1 = 1,..., 4, 



(83) 



fc=i 



where 



(84) 



ATW = 2xsin(fx)T [(A||l«)/^(a,) + (E\\E®)/S A {u) + (f\\T®)/S t (w) 

and where M"" k ' is the 4x4 matrix with components 

M (0(fe) = 4x 2 sin 2 (§x)T \(A^\\A^)/S a (lu) + (E^\\E^)/S a (lu) + {f^\\f^)/S f (uj)\ . (85) 

The solution of Eq. (|83")) is simplified by noticing that the component functions A^(t), E^ k \t), and T^ k \t) consist 
of simple sines and cosines with period ~ 2-7r/w, modulated by the slowly changing functions Ui(t) and Vi(t) (with 
periods that are multiples of 1 yr). By the approximate orthogonality of sine and cosine terms, for To ~ 1 yr the 
scalar products (A^ \ \A">) can be approximated as 



and 



(AW\\aW)~(AP>\\AW)~0, 



(AWpW) ~ (A^\\A^) = \U A , 
(A^\\A^) ~ {A^\\A^) = \V A , 



(86) 



(87) 



1 P- 



(l«pW) ~ -(A( 2 )p( 3 )) 

with similar expressions for the E^ and T^ k \ The matrix then simplifies to 



(89) 
(90) 



M (0W = Z£ 



( U Q P\ 

Q V -P 

-P [7 Q 

\P Q V J 



(91) 
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where the elements U, V , Q, and P are given by 



U = 4x 2 sm 2 a X )[U A /S A (uj) + U E /S A (uj) + U T /Sf(u>)} 



(92) 

V = Ax 2 sin 2 (|x) [V A /S A {u) + V e /S a (lj) + V?/Sr(u>)] , (93) 
Q = 4x 2 sin 2 (|ar) [Q A /S A (lj) + Q E /S A [u) + Qf/S^)] , (94) 
P = 4x 2 sin^fz) [P a /S a (uj) + P E /S A (u) + Pt/S t {uj)] . (95) 

Then the analytic expressions for the maximum likelihood estimators of the amplitudes a^ k ' are given by 



(96) 



where A = UV - Q 2 - P 2 . 

Substituting the ML amplitude estimators a^' in the likelihood function A yields the reduced likelihood function 
A r . The logarithm of A r is known as the T statistic; using Eqs. (|79[). ([8lTl . ([84]) . and (|96|) . we find 









V 


-Q 









( \ 




2 




-Q 


u 


p 







JV(2) 











p 


r 


-Q 




JV(3) 


v « (4) y 






-p 





-0 


u J 




V A^ 4) / 



^ = 3 EE (M-^^^ATWiVW = (ToA)- 1 {y [(iV«) 2 + (M 3 >) 2 ] + U [(N^) 2 + 



=i fe=i 



2Q [jV^W 2 ) + iV^ 4 )] - 2P [n^N^ - N m N (3)'] I _ (97) 



We adopt T as the detection statistic of our proposed search scheme. The statistic T is already maximized over 
the amplitudes cS k \ which are known in this context as extrinsic parameters. By contrast, the ML estimators of the 
remaining [intrinsic) source parameters (u>, u>, (3, and A) are found by maximizing T . In practice, this is done by 
correlating the detector output with a bank of optimal filters precomputcd for many values of the intrinsic parameters. 
Introducing the complex quantities 

(98) 
(99) 
(100) 
(101) 
(102) 



W 

N (v) 



,(2) 



in 



(3) 

) 

,(4) 



Q + iP, 



we can write the ML amplitude estimators and the T statistic in the compact form 



[where A = UV - \W\ 2 } and 

T 



a (u) = 2(T A)- 1 [WV (u) - W*NW], 
a {v) = 2{T A)- 1 [UN^ -WN^], 

(T A)- 1 {t/|A^ (tl) | 2 + C/|iVW| 2 -2Re WN^(N^)* } 



(103) 
(104) 

(105) 



In Sec. |V] we shall see that this expression is very suitable for numerical implementation. Equations (| 1 03j) (|105|) 
summarize the proposed ML data-analysis scheme, which uses all the available LISA data. Similar expressions 
hold if we analyze a single interferometric combination, such as X\. In Appendix fAl we describe a useful complex 
representation of the GW TDI responses that simplifies the integrals involved in the computation of T and of the ML 
amplitude estimators. 

In the LW approximation,_Eqs. ([105 ]) — ([105 ]) simplify somewhat: using the A^ (t) and E$$(t) of Eqs. (JBDJ and 
([6Tj) [and remembering that T^ff(t) ~ 0], we go through with our formalism in parallel with Eqs. (|84l) (f9"Tj) . and find 
that Plw — 0, so Wlw ^ s rea l- The complex variables and N^ v ' are given by the integrals 



Ntt = -2i {UJL) 



""SW( U )J D 



N& = -2i {UJL) 



^V2[2u 2 (t) - m(t) - u 3 (t)]A(t) + 3V6[u 3 (t) - ui(t)]E(t)\ e #(t) dt, (106) 
{3V2[2w 2 (i) - Vl (t) - v 3 {t)]A(t) + 3V6[v 3 (t) - Vl (t)]E(t)X e**® dt. (107) 



Analogous LW expressions hold for a single TDI observable, such as X±. 
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B. Distribution of the T statistic 



Crucial to a search scheme based on comparing the ML statistic T with a predefined threshold is the determination 
of the false-alarm probability Pp (which determines how often T will exceed the threshold in presence of noise alone) 
and of the detection probability Pjj (which determines how often T will exceed the threshold when a signal is present, 
resulting in correct detection). In this section we compute the probabilities Pp and Pp> for the correlation of detector 
data against a single optimal filter (i.e., for fixed values of the intrinsic parameters). 

Under the assumption of zero-mean Gaussian noise, the weighted correlations [Eq. (|84|) ] are Gaussian random 
variables; since T is a quadratic form in the [see Eq. (1^71) ]. it must follow the % 2 distribution. Following Sec. 

Ill B of Ref. [22j], we can diagonalize the quadratic form to find that, in the absence of the signal, IT follows the x 2 
distribution with n — 4 degrees of freedom [29(. In presence of the signal, IT follows a noncentral x 2 distribution 
with n = 4 degrees of freedom and with noncentrality parameter k equal to the optimal (S/N) 2 = p 2 [see Eq. ([76]) ]. 
For instance, if we use A, E, T, 

k = p 2 = (A GW \A GW ) + (E GW \E GW ) + (T GW \T GW ) (108) 

(which agrees with the result derived in Ref. [2l| ) , while if we use only Xi , 

K . = p 2 = (X GW \X GW ). (109) 



The x 2 probability density function is 



j:n/2-l 

P0(f)= (^ HPK) (110) 



for k = 0, or 



Pl(p;F) = ^n'n-l"' In/2-l(pV0) CX P (-^- \p 2 ) (111) 



(n/2-l)/2 



for k = p 2 , where I n /2-i is the (n/2 — l)th-order modified Bessel function of the first kind. Thus, the false-alarm 
probability for a threshold JF is 

n/2-l 



P F {F )= / p (^)^=exp(-jr ) Y, rf/kl (H2) 



(for even n; for odd n the result involves the error function) while the detection probability, in the presence of a 
S/N = p signal (using the correct optimal filter), is 



/>oo 

Pd(p;Fo)= / pi(p,J r )dJ : ; 



(113) 



this integral cannot be evaluated in closed form in terms of known special functions, but it is clear that the higher 
the optimal S/N, the higher the detection probability. 



C. False-alarm and detection probabilities for GW searches 



In actual GW searches, the detector output will be correlated to a bank of optimal filters corresponding to different 
values of the intrinsic parameters £ M . For a given set of detector data, the statistic ^-"(^ M ) is a generalized multipa- 
rameter random process known as random field (see Adler's monograph 30] for a comprehensive discussion): we can 
use the theory of random fields to get a handle on the total false-alarm and detection probabilities for the entire filter 
bank. 

We define the autocovariance C of the random field ^ r (^ M ) as 

c(£ M ,r) = £o[He)H^)}-£o[He)}£o[H^)}, (114) 

where the expectation value £q is computed over an ensemble of realizations of noise (in absence of the signal). In Ref. 
[IH the total false-alarm probability was estimated by noticing that the autocovariance function tends to zero as the 
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displacement A£ M = — £ M increases (and in fact, it is maximum for A£ M = 0). The space of intrinsic parameters 
may then be partitioned into a set of elementary cells, whereby the autocovariance is appreciably different from zero 
for within each cell, but negligible between cells. The number of elementary cells needed to cover the parameter 
space gives an estimate of the number of independent realizations of the random field (i.e., the number of statically 
independent ways that pure noise can be strongly correlated with one or more of the optimal filters) . 

There is of course some arbitrariness in choosing the boundary of the elementary cells; wc define them by requiring 
that the autocovariance between the center and the surface be one half of the autocovariance at the center: 

C(f, = ±C(£", £ M ) for at cell center, £ /Al on cell boundary. (115) 

Taylor-expanding the autocovariance to second order in A£^, we obtain the approximate condition 



2 <9£' p ar 



A^Ar = §C(£",£"), (116) 



with implicit summation over p and a. Within the approximation (necessary to obtain results in simple analytical 
form), the cell boundary is the (hyper-)ellipse defined by G pCT A£ /> A£ <7 = 1/2, where [3l[ 



i i o 2 c{^,e) 



(117) 



(in Appendix |5] we shall derive a relation between this G pa and the Fisher information matrix) . The volume V cc \\ of 
the elementary cell is then 

(7T/2)*/' 

Fcc " " r(^/2 + l)VdetG p / (118) 

where if is the number of intrinsic parameters, and T is the gamma function. The total number of elementary cells 
within the parameter volume V is given by 

Neon = r ( ( f / y / 2 ) J V^GpZdV. (119) 

As discussed above, we consider the values of the statistic T within each cell as independent random variables, which 
in the absence of signal are distributed according to Eq. (|110p . By our definition of false alarms, the probability 
that T will not exceed the threshold T in a given cell is just 1 — Pf(Tq); the probability that T will not exceed the 
threshold Tq in any of the cells is 

1 - PF.tot^o) - [1 - Pf(Fo)] Nm ] (120) 

this Pp, tot (-^o) is therefore the total false- alarm probability for our detection scheme. 

When the signal is present, a precise calculation of the probability distribution function of T is nontrivial, since 
the presence of the signal makes the random process x(i) nonstationary. However, we can still use the detection 
probability given by Eq. (|113[) for known intrinsic parameters as a substitute for the detection probability when the 
parameters are unknown. This is correct if we assume that, when the signal is present, the true values of the intrinsic 
parameters fall within the cell where T is maximum. This approximation is accurate for sufficiently large S/N. 

V. FAST COMPUTATION OF THE T STATISTIC 

The detection statistic T [Eq. (| 105[) ] involves integrals of the general form 

x(t) m(t; u>, P, A) exp[i<fi mo d(t; to, to, (5, A)] exp[iw£] dt (121) 

o 

where m is a combination of the complex modulation functions defined in App. A, while the phase modulation ^> mo d 
is given by 

4>mod{t) = T^tot 2 + ooR cos fi cos(fM + ?] ~ A) (122) 
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[Eq. (|45p]. We see that the integral (|121[) can be interpreted as a Fourier transform [and computed efficiently with a 
fast Fourier transform (FFT)], if <j) mo d and m do not depend on the frequency u>. In fact, even in that case we can 
still use FFTs by means of the procedure that we now present. 

From the original data we generate several band-passed data sets, choosing the bandwidth of each set so that 
m exp[i(f) mo d\ is approximately constant over the band. We then search for GW signals in each band-passed data set: 
this is done by computing the T statistic over a grid in the parameter space (w, (3, A), set finely enough that we do not 
miss any signal. We follow the grid-construction procedure presented in Sec. Ill A of Ref. (32J. The phase modulation 
can be usefully reparametrized as 

4>mod = Pit 2 + Acos(ftt) + Bsm(Q,t), (123) 

where 

Pi = 

< A = uR cos (3 cos(\- rj ), (124) 
B = cji?cos/3sin(A — 770). 

Since m is a slowly changing function of time, we consider it constant for the purpose of constructing a grid over 
the parameter space. The result is a uniform grid of prisms with hexagonal bases, where the parameter subspace 
A-B is tiled by regular hexagons. The grid in the parameters Co, (3 and 7 is then derived by applying the inverse 
transformation, 

!Co = 2pi, 
(3 = ± arccos [^A 2 + B 2 /uR s j , (125) 
A = rjct + arctan(£?/A), 

where for each band-passed data set we set the unknown frequency lo to the maximum frequency of the band. The 
computation of the T statistic includes both phase- and amplitude-modulation effects, even if these were neglected 
in the construction of the grid [in fact, the sign degeneracy for (3 in Eq. (|125p is resolved by amplitude modulation, 
which distinguishes between sources in opposite directions with respect to the plane of the ecliptic]. 

Once we have a detection, the accurate estimation of signal parameters requires a second step. Since the coarse 
signal search described above is performed by evaluating the function m exp[i(j) mo d\ at the maximum frequency of each 
band, our filters are not perfectly matched to the signal, and thus are not optimal; as a consequence, the location 
of the maximum of T does not correspond to the correct ML estimators. We therefore refine the coarse search by 
maximizing T near the coarse-search maximum, this time without any approximation. 

We have performed a few numerical simulations to assess the performance of our optimal-filtering algorithm. Here 
we report on three of them. In the first simulation we analyzed the X\ TDI data corresponding to two simultaneous 
monochromatic signals, of frequency / = 3 mHz and S/Ns of 24 and 10, emitted from sources at opposite positions 
with respect to the plane of the ecliptic. We generated a one-year-long time series for X\ by implementing Eq. (|58[) 
numerically and we included noise by adding a Gaussian random process (as realized by a random number generator) 
with spectral density given by Eq. (|65p. We narrowbanded the X\ data to a bandwidth of 0.125 mHz around 3 mHz, 
and we analyzed the resulting data by implementing the two-step procedure described above, using the Nelder-Mead 
maximization algorithm [33| for the second step. The angular grid for the all-sky search consisted of about 900 
points. We then performed the following operations: i) detecting the stronger signal and estimating its parameters; 
ii) reconstructing the stronger signal and subtracting it from the data; iii) detecting the weaker signal and estimating 
its parameters; iv) subtracting it from the data. Figure [5] shows the amplitude spectrum of X\ before and after the 
subtraction of the two signals, as compared with the spectrum of noise alone. Figure [7] shows a comparison of the 
input signals with the reconstructed signals (built with the parameters specified by the ML estimators). We see that 
the amplitude modulations in the GW response enable us to determine the sky location of two sources of the same 
frequency, and also to resolve the two GW polarizations. Signal resolution will degenerate rapidly as more sources of 
the same frequency are added, so the steps described above cannot be used as a general signal subtraction procedure 

In the second simulation we analyzed the X\ TDI data corresponding to a single signal of frequency / = 25mHz, 
S/N = 9.5, and / = 6.5 x 10~ 13 Hz s _1 (corresponding to a binary of chirp mass M c = O.9Af ). We generated a 
one-year-long time series for X\ by implementing numerically the exact GW response, Eq. (|48[) . and we added noise 
as described above. We narrowbanded the X\ data to a bandwidth of 0.5 mHz, and again we analyzed the resulting 
data with the two-step procedure described above. The sky search was performed on a small grid (~ 300 gridpoints) 
around the true values of the signal parameters. In the third simulation we analyzed the A, E, and T TDI data 
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FIG. 6: Quality of signal reconstruction, as seen in the Fourier domain, in the first simulation. The original X\ time series 
contains noise, pluse two monochromatic signals of equal frequency (3 mHz) and ecliptic longitude, but opposite ecliptic 
latitudes. We show the X\ amplitude spectrum before and after the subtraction of the reconstructed signals, compared with 
the amplitude spectrum of noise alone. 



corresponding to the same signal, for a total S/N = 19. The ML search procedure was performed as in the second 
simulation. The top panel of Fig. [5] shows the TDI observable A(t) for the signal alone, superimposed on the TDI 
observable for signal plus noise: we see that the signal is more than one order of magnitude weaker than the noise. 
The bottom panel of Fig.[S]shows the T statistic (already maximized over /, f3, and A) near the input signal frequency. 
We see that the statistical significance is higher for the multiple-observable search than for the X\ search. Figured 
shows a comparison of the input signals with the reconstructed signals. We see that reconstruction is more accurate 
for the multiple-observable search, but in both searches our procedure resolves the two GW polarizations successfully. 

We conclude that our proposed algorithm performs satifactorily, detecting the simulated signals, accurately esti- 
mating their parameters, and resolving the two GW polarizations, both in the low-frequency regime (first simulation) 
and high-frequency regime (second and third simulation). In a future paper, we plan to discuss in detail the expected 
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FIG. 7: Quality of signal reconstruction, as seen in the time domain, in the first simulation. The panels show the input signals 
(the stronger on the left, the weaker on the right), compared with the reconstructed signals; the two GW polarizations are 
plotted together (top row), and separately (middle and bottom rows). 



errors in parameter estimation for a source with given frequency, sky position, and S/N. 
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APPENDIX A: COMPLEX REPRESENTATION OF THE RESPONSE 

From Eqs. (|52|) and |53|) it is easy to see that the Sagnac TDI observables ati can be rewritten in the complex form 



af w = 2a;sin(§z) Re [a (u) *m£f (t) e**<*> + a {v) *m^(t) e**<*> 



(Al) 



where x = ljL, the complex amplitudes and are defined in Eq. I|98p and the phase <f>(t) by Eq. (jlS)) . and 
the complex modulation functions rtia) and are given by 



m { Z ] (t) 



E 

3=1 



Uj(t) 
Vj (t) 



e- txd > {sinc[(l + Cj )x/2]n+i + sinc[(l - Cj )x/2]n^} 



(A2) 
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FIG. 9: Quality of signal reconstruction, as seen in the time domain, in the second and third simulations. The panels show 
the input signals (Xi on the left, A on the right), compared with the reconstructed signals; the two GW polarizations are 
plotted together (top row), and separately (middle and bottom rows). Signal reconstruction is more successful for the multiple- 
observable search (right) than for the single-observable search (left). 



with Uj(t), vj(t) given by Eqs. (f27|) and ([28]) . Cj(t) and dj(t) by Eqs. (|47|) and (|46]) . and with the constants given 
in the left part of Table HO 
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TABLE II: Constants that appear in the complex representation of the GW responses of the TDI observables. The constants 



71^2 and n^l are obtained from n^l by cyclical permutation of the index j, as are n^ 2 and from nj^. 

The optimal combinations A GW , E GW , T GW are given by formulas similar to Eq. (|A1[) . with modulation functions 



™(«) 

A 



(u) 
m E 



V2 V w " 3 

— fm (u) 



,(«) 



'3 J ' 



(A3) 
(A4) 
(A5) 



and similar expressions for m ( j\ m { ^\ and . The quantities N<- u \ , U, V, and W [Eqs. (pUT]) . (fItJ2"j) . 

(jlOOp ] that are needed to compute the ML amplitude estimators a"' and the statistic [Eq. (| 105|) ] can be 
written in terms of the complex modulation functions as 



A(t) mf(t) + E(t) mf(t) T{t) m£°(t) 



and 



iV W = 4a;sin(fa;) 



U = 4a; 2 sin 2 (fa;) (2/T ) 



V = 4a; 2 sin 2 (§a;)(2/To) 



W = 4a; 2 sin 2 (fx) (2/T ) 



'ih 



A(t) m%Xt) + E{t) m { ;\t) T(t) m { ;\t) 



S a {uj) 



e^dt, 



|m ( j%)| 2 



S A {u) 
\m^\t)\ 2 + \m^W 



Sf(u>) 



d.t, 



dt, 



Nto = 4a;sin(|a;) / " v "' "" A ^ - ^ "" E w + : ""u" e^dt, (A6) 

(A7) 

(A8) 
(A9) 

dt. (A10) 
(All) 

the modulation functions (t) and m x (t) have exactly the same functional form as the functions mi"\f), rria}(t) 

defined in Eq. (|A2[) , except that the coefficients are those given in the right part of Table [TTJ For the single Xi 
observable, the ML estimators for the amplitudes and for T are again given by 



S a {uj) S t (uj) 
mf\t)m ( l\t) + mf\t)mf{t) m^\t)mf{t) 



S a {uj) 



The Xi TDI observables can be written in the complex form as 

X GW = 4a;sin(x) sm(2x) Re \ia {u) * m [ ^{t) e^ (t) + ia {v) * m%]{t) e^ (t) 



a« = 

a(<0 - 



2(T A Xl )-i[V Xl NW -W* Xl N$; 



2(T A Xl )-i[U Xl N<£-W Xl NW] 



r(«)l 



[where A Xl = E^V^ - l^xj 2 ] and 

^ = (TpA)- 1 {v^ lA^I 2 + U Xl \N%\ 2 - 2 Re N^(N%)*] } , 
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APPENDIX B: REDUCED INFORMATION MATRIX 



It is interesting to examine the relation between the matrix G^ v defined by Eq. (|117[) and the Fisher information 
matrix T %3 . We consider the case of a single TDI observable; multiple observables can be treated in similar fashion. 
As seen in Sec. II, the generic TDI GW response h(t) can be written as the linear combination 



(Bl) 



k=l 



as discussed in Sec. IIV A[ the amplitudes are extrinsic parameters, while all the other parameters (denoted 
together as £ M ) are intrinsic [all the parameters are denoted together as 6 l = {a^ k \^)]. Note that in the case 
of the TDI observables Ap w or ^ GW [Eqs. l[3g|). l[8"T])]. the component functions h^(t) would include the factors 
Axsm(x) sin(2x) and 2xsin(|a;) respectively. 

In this notation, it is easy to show that the optimal S/N and the Fisher matrix can be written as 



a T ■ M • a, 



and 



M 



Fa 



Sa 



(B2) 



(B3) 



where the top and left blocks correspond to the extrinsic parameters, while the bottom and right blocks correspond to 
the intrinsic parameters. The superscript " T " denotes transposition over the extrinsic parameter indices. Furthermore, 
a = (a' 1 ', a,( 2 \ a^ 3 \ a*- 4 -*), and the matrices M, F, and S are given by 



MWW = (fe<*>|fcW) 



S iW) = (^H 



(B4) 
(B5) 

(B6) 



The covariance matrix C*- 7 , which expresses the expected variance of the ML parameter estimators, is defined as 
(r _1 ) lJ . Using the standard formula for the inverse of a block matrix [35[ we have 



M- 1 + M- 1 ■ (F ■ a) ■ f- 1 • (F • a) T • M" 1 
-f- 1 • (F-a) T • M- 1 



- 1 • (F-a) -f- 1 
f-i 



(B7) 
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where 

f = a T • (S- F T • M" 1 • F) • a. (B8) 

We shall call (the Schur complement of M) the projected Fisher matrix (onto the space of intrinsic parameters). 
Because the projected Fisher matrix is the inverse of the intrinsic-parameter submatrix of the covariance matrix C tJ , 
it expresses the information available about the intrinsic parameters once the extrinsic parameters are set to their 
ML estimators. Note that is still a function of the putative extrinsic parameters. Using Eq. (|B2|) we define the 
normalized projected Fisher matrix 

r * sI V = a r. M . a (B9) 

From the Rayleigh principle [35| . it follows that the minimum value of the component F^" is given by the smallest 
eigenvalue (taken with respect to the extrinsic parameters) of the matrix [(S — F T • M" 1 • F) • M _1 ] MI/ . Similarly, the 
maximum value of the component T^ 1 is given by the largest eigenvalue of that matrix. Because the trace of a matrix 
is equal to the sum of its eigenvalues, the matrix 

f = ^-Tr [(S — F T • M _1 ■ F) • M _1 ], (BIO) 

where the trace is taken over the extrinsic-parameter indices, expresses the information available about the intrinsic 
parameters, averaged over the possible values of the extrinsic parameters. Note that the factor \ is specific to the case 
of four extrinsic parameters. We shall call the reduced Fisher matrix. This matrix is a function of the intrinsic 
parameters alone. 

Let us now compute the components of G^, defined by Eq. (|117[) . We start from log A, which in our notation is 
given by 

logA = a T -N-±a T -M-a, (Bll) 

where — (x\h^ k '), with x(t) — n(t) + h(t), and n(t) a zero-mean Gaussian random process. The ML estimators 
a^ are given by a = M _1 • N, so for the T statistic we have T — \ N T • M _1 • N. Using the relations [13] 

£{(n|ai)H* 2 )} = (*i|*2), (B12) 
E {(n|si)(n|s 2 )(n|s 3 )(n|s4)} = (si|s 2 )(s 3 |s 4 ) + (si|s 3 )(s 2 |s 4 ) + (si|s 4 )(s 2 |s 3 ), (B13) 

where s\, s 2 , s 3 , and S4 are deterministic functions, we find that the autocovariance function C(^,^ /AI ) of Eq. (|114p 
is given by 

C(e, £' M ) = ^Tr [Q T • M- 1 • Q • M' -1 ], (B14) 

where 

Qui = {hW\tiW), (B15) 

and the primes denote functions of the primed parameters Inserting Eq. (|B14|) into Eq. (jTT7J), after some lengthy 
algebra (omitted here) we come to the final result 

G ^ = \ E [S« (i) -i^ m)(fe HM- 1 )( m )Wi^")W](M- 1 )WW =IV (B16) 

k,l,m,n 

Thus, the JF-statistic metric G A " / [3l| is found to be exactly equal to the reduced Fisher matrix f that this should 
be the case is understandable, since both matrices contain information about the relatedness of waveforms with 
nearby values of their intrinsic parameters (while both assume that the extrinsic parameters arc being set to their ML 
estimators). For a related argument about the placement of templates for a partially maximized detection statistic, 
see Ref. 1361. 
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APPENDIX C: FIRST-GENERATION TDI RESPONSES 

The GW response of the first-generation TDI observable X is given by Q 

yGW _ C,,GW , „GW\ , /.,GW , ,,GW\ /..GW , „GW\ (..GW , ..GW) . 

A - 12/31 +2/13,2) + (2/21 + J/l2,3 J, 22 - (2/21 + 2/12,3) - 12/31 + 2/13,2 ),33, 

after some algebra we get to 

4 
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«a(t) 






_«2(t)_ 



X GW (4) = sin(wL)^ a Wx( fe )(/), 
fe=i 

where the functions X( fe '(4) are given by 

{sinc[(l + c 2 )x/2] cos [0(f) — xd,2 — 3x/2] + sinc[(l — c 2 )x/2] cos [0(4) — xd2 ~ 5x/2] } 

{sinc[(l + c 3 )x/2] cos[0(4) - xd 3 - 5x/2] + sinc[(l - c 3 )x/2] cos [0(4) - xd 3 - 3x/2] }(C3) 

{sinc[(l + c 2 )x/2] sin [0(4) - xd 2 - 3x/2] + sinc[(l - c 2 )x/2] sin [0(4) - xd 2 ~ 5x/2] } 

{sinc[(l + c 3 )x/2] sin [0(4) - xd 3 - 5x/2] + sinc[(l - c 3 )x/2] sin [0(4) - xd 3 - 3x/2] },(C4) 



«3(*) 



" X (3)~ 




W 2 (4) 






_v 2 (4)_ 



«s(t) 
«s(*) 



where x = wL and sine . . . = (sin ...)/(...). The GW responses for Y and Z can be obtained by cyclical permutation 
of the spacecraft indices. 

The GW response for a can be written in similar form: 



fe=4 



a GW (t)=u;LY,a {k) u {k) (t), 



(C5) 



fc=i 



where 







«i(t) 


a( 2 ) 







a 



«( 3 ) 
(4) 



tta(t) 
«a(t) 

u 3 (t) 

«3(*) 

Ul(t) 
«l(t) 

ua(t) 
«a(*) 

«a(t) 



{sinc[(l + Ci)x/2] sin [0(4) - xd x - 3x/2] - sinc[(l - c x )x/2] sin [0(4) - xd\ ~ 3x/2] } 
{sinc[(l + c 2 )x/2] sin [0(4) - xd 2 - x/2] - sinc[(l - c 2 )x/2] sin [0(4) - xd 2 - 5x/2] } (C6) 
{sinc[(l + c 3 )x/2] sin [0(4) — xd 3 — 5x/2] — sinc[(l — C3)x/2] sin [0(4) — xd 3 — x/2] }, 
{sinc[(l + ci)x/2] cos[0(4) — xd\ — 3x/2] — sinc[(l — ci)x/2] cos[0(4) — xd\ — 3x/2] } 
{sinc[(l + c 2 )x/2] cos [0(4) - xd 2 - x/2] - sinc[(l - 03) x/2] cos [0(4) - xd 2 - 5x/2] } (C7) 
{sinc[(l + 03)21/2] cos [0(4) — xd 3 — 5x/2] — sinc[(l — C3)x/2] cos [0(4) — xd 3 — x/2] }. 



The and 7 combinations are again obtained by cyclical permutation of the spacecraft indices. 
For the TDI observable £ we find 



C GW (4) = cL^ a W C W(4), 



(C8) 



fc=i 
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with 



- C (l>- 




Ui(t) 


((2) 




Vl (t)_ 



+ 



{sinc[(l + c\)x/2\ sin [0(f) — xd% — 3x/2] — sinc[(l — ci)x/2] sin [0(f) — xd\ — 3x/2]} 



u 2 (f) 
v 2 (t) 

u 3 (t) 

C (3) l = _ k(t) 

C (4) wi(t) 
u 2 (f) 

W2(<) 
«3(*) 

«s(t) 

Finally, the optimal TDI observables A, E, and T [2l| are defined as linear combinations of a, (3 and 7: 

A = ±fr-a), 

£ = -^(a-2/3 + 7), 

T = _L(a + /3 + 7 ). 



{sinc[(l - 


V c 2 )x/2] 


sin [0(f) 


- xd 2 


- 3x/2] 


— sinc[(l 


- c 2 )x/2] 


sin [0(f) - 


- xd 2 


{sinc[(l - 


\- c 3 )x/2] 


sin [0(f) 


- xd 3 


- 3x/2] 


— sinc[(l 


- c 3 )x/2] 


sin [0(f) - 


- xd 3 


{sinc[(l - 


hci)i/2] 


COS [0(f) 


— xd\ 


- 3x/2] 


— sinc[(l 


-a)x/2] 


COS [0(f) 


— xd\ 


{sinc[(l - 


Vc 2 )x/2] 


COS [0(f) 


- xd 2 


- 3x/2] 


— sinc[(l 


- c 2 )x/2] 


COS [0(f) 


- xd 2 


{sinc[(l - 


\- c 3 )x/2] 


COS [0(f) 


- xd 3 


- 3i/2] 


— sinc[(l 


- c 3 )x/2] 


COS [0(f) 


- xd 3 



(Cll) 



The long-wavelength approximation to the GW responses is obtained by taking the leading-order terms of the 
generic expressions in the limit of ojL — > 0. For instance, for X [Eqs. (|C2[) - (|C4[) ], we get 



X 



GW 
LW 



A{loL) 2 |[u 2 (f) - u 3 (t)] a (1) cos 0(f) + a (3) sin 0(f) + [v 2 (t) - v 3 (t)\ a (2) cos 0(f) + a (4) sin 0(f) } , (C12) 



with given by Eqs. (|3T ]) -(|3I |) . and m(t), v % {t) by Eqs. (J27J) , (J2HJ) • The LW responses for Y and Z can be obtained 
by cyclical permutation of the indices. Adopting the notation of Ref. Q, we find also that 



X™(t) * 2L 2 (n 3 L )' • H L (t) • - (n 2 L )' • H L (f) • 



(C13) 



where """ denotes the second time derivative, is given by Eq. and H L (f) is given by Eq. (fill) . 

The GW response of the Sagnac observable a^w is equal simply to ^X^ . From Eqs. (jClip we then get the LW 
GW responses A§™, Eg$, and T§™: 



iGW 



V2(luL) 2 U-2u 2 (t) +«i(f) +u 3 (t)] a (1) cos0(f) +a (3) sin0(f) 



+ [-2v 2 (t)+Vi(t)+v 3 (t)] a (2) cos 0(f) +a (4) sin 0(f) |, 



(C14) 



V6(loL) 2 |[in(f) -u 3 (*)] a (1) cos 0(f) +a (3) sin 0(f) + [ui (f) - « 3 (f )] a (2) cos 0(f) + a (4) sin 0(f) }, (C15) 



0[(luL) 3 



(C16) 
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